#include <iostream>
#include "Okada.hxx"

int main()
{
  int n=100;
  double size=10;
  double dx=size/n;

  double lambda(1.2), mu(2.3), slip(-10), opening(0), width(2.5), length(2),
    strike(0), dip(90), rake(0);
  FTensor::Tensor1<double,3> origin;
  origin(0)=size/2;
  origin(1)=size/2;
  origin(2)=6;

  Okada okada(lambda, mu, slip, opening, width, length, strike, dip,
              rake, origin);

  for(int i=0;i<=n;++i)
    for(int j=0;j<=n;++j)
      {
        FTensor::Tensor1<double,3> coord(i*dx,j*dx, size/2);
        Okada::Displacement d(okada.displacement(coord));
        std::cout.precision(15);
        std::cout << std::scientific
                  << coord(0) << " "
                  << coord(1) << " "
                  << coord(2) << " "
                  << d.first(0) << " "
                  << d.first(1) << " "
                  << d.first(2) << " "
                  << d.second(0,0) << " "
                  << d.second(0,1) << " "
                  << d.second(0,2) << " "
                  << d.second(1,0) << " "
                  << d.second(1,1) << " "
                  << d.second(1,2) << " "
                  << d.second(2,0) << " "
                  << d.second(2,1) << " "
                  << d.second(2,2) << " "
                  << "\n";
      }
}
